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The response of the ocean to an axisymmetric cyclone is 
investigated using a model based on the conservation of potential 
vorticity in the interior. The results, which confirm the conservation 
potential vorticity, are used to obtain limiting values of the 
vertical eddy diffusivity and an effective Ekman depth. Some aspects 
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I. INTRODUCTION 


In October of 1968, Tropical Storm Kara forced R/V PILLSBURY 
south from an operating area in the vicinity of latitude 30°925'N, 
longitude 76°52'W. While in the stages of maturation, the tropical 
storm displayed the curious behavior of slowing and crossing its own 
path several times. Kara eventually reached hurricane force, increased 
speed and moved on a generally northeast course. In transit to the 
Operating area after passage of the storm, the research vessel 
experienced set to the right and then left, indicating that a hurricane- 
induced eddy had been crossed. Bathythermograph stations were taken 
in order to locate the associated dome-like density structure,and STD 
stations across the presumed dome were taken after it had been located. 
The track of Kara and the track of R/V PILLSBURY are shown in Figure l. 
The hurricane track was obtained from the advisories issued by the 
National Hurricane Center. 

Kara may be considered to have been quasi-stationary during the 
period 14-15 October 1968. Therefore it affords an opportunity to 
study some aspects of hurricane-induced dynamics through the upper 
part of the water column. As Stevenson (1966) points out, hurricanes 
are unique "laboratories" for the study of air-sea interaction and 
ocean dynamics. Although the significant atmospheric effects of a 
hurricane are spread over horizontal scales on the order of several 
hundred kilometers, the most interesting oceanic effects take place 


within an area which can be crossed in one day. Unfortunately, 
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investigations of the response to hurricanes depend by necessity on 
"cruises of opportunity." Given the availability of research vessels, 
studies of this sort are further complicated by the unpredictability 
of the occurrence, duration and track of a storm. 

A study of low intensity storms such as Kara is particularly 
relevant because the wind velocities are not too different from those 
found in the winter storms in the North Atlantic Ocean. In particular, 
the data provide a means of testing assumptions of conservation of 
potential vorticity in the ocean. An analysis of those areas where 
the conservation equations are not applicable is important because 
understanding the causes of the differences between observed and 
measured response contributes to an understanding of the response of 
the ocean in the near surface layer. A third reason for studying such 
storms is that the deep ocean response may be typical of numerous 
storms of near hurricane force. If so, the associated changes in mass 
structure are significant at great depths. For example, the 
perturbation of the isopycnals will have drastic effects on acoustic 
propagation at depths at least as great as the SOFAR channel. The 
analysis also serves as an addition to a catalogue of techniques useful 
in analyzing transient oceanic behavior. 

Circular symmetry is assumed in the model used in this paper. 
However, it can be seen from Figure 1 that insufficient stations have 
been taken to positively identify the center of the vortex. Station 24 
was made after the vortex had been crossed once, when the research 
vessel returned to take a center station. It was felt that the station 


was taken close to the center but the possibility that it was not 





cannot be excluded. The stationarity of Kara increases the amplitude 
of the response. Since the data were collected several days after the 
storm passed,it may be assumed that most of the energy is found in 


geostrophic motion. 





IT. METHODOLOGY 


The existence of a cold water dome after the passage of a hurricane 
is suggestive of upwelling. This was concluded by Leipper (1967) and 
is supported by data collected by Landis (1966). Leipper's data 
inspired O'Brien and Reid (1967) to pursue their study of upwelling in 
a two layer ocean as a response to a stationary hurricane. O'Brien 
and Reid established through a transient, non-linear model that the 
wind stress associated with storms of hurricane force and structure is 
capable of producing velocity divergence in the surface layer and large 
scale upwelling. O'Brien (1967) added turbulent mixing to the model 
and produced results which bear resemblance to Leipper's observations 
of the near surface layer. 

The problem can be approached somewhat differently. It is well 
known that surface divergence and associated upwelling will stretch 
water columns. A fluid will respond by spinning, thereby generating 
relative vorticity about the vertical axis. The effect extends to 
the bottom of the fluid and is transmitted much more efficiently than 
momentum transfer by friction. The model in this study utilizes the 
conservation of potential vorticity to describe the interior response. 
The assumption of a frictionless interior makes the continuously 
stratified problem tractable by allowing use of the conservation 
equation in the interior. This model and O'Brien's differ because 
they investigate different aspects of the response to axisymmetric 
forcing. O'Brien solves the governing equations numerically as an 
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initial value problem in time because he is interested in the evolution 
of the response. By its very nature, his two layer model is incapable 

of obtaining the structure of the interior. However, it is successful 

in duplicating many aspects of the time dependent response. This 

study seeks to model the steady-state response of the ocean independently 
of the generating mechanism using a realistic density profile. 

The model used in this study was derived by Rooth (1970) and is 
presented fully in Appendix B. It makes use of the non-linear, time 
dependent equations of motion in deriving the vorticity equation. 
Hydrostatic equilibrium is assumed for the formation of the divergence 
equation. It is assumed that the data were collected a sufficient 
length of time after the storm passed to ensure that motion was 
essentially geostrophic. Combining the divergence equation and the 
vorticity equation results in a governing equation expressing the 


conservation of potential vorticity in a geostrophic velocity field, 


£2 92e 


y2E + = 0. =| 
: N2(Zp) %Zo2 


€ is the lifting of a density surface from its undisturbed state (called 
the matching level), f is the Coriolis parameter, and N2(Zo) is the 
square of the Brunt-Vaisala frequency of the undisturbed state. Zg is 
the geometric height in the undisturbed state. The observed mass 
distribution perturbed from the undisturbed state. The equation is 

the result of a transformation from cartesian coordinates with geometric 
height as the independent variable in the vertical to a "natural" 
coordinate system where density is the vertical coordinate. Figure 2 
illustrates the relationship between § and Zp. Geometric depth is the 


sum of Zp and §. Formally, Z = Zo(p) + § (Zo.r), where r is the radius. 
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Figure 2. The relation between density (?), undisturbed depth (Zo) and 
Siretching (&). 





Not only is the numerical solution considerably simplified in this 
formulation, but solving in terms of the stretching is intuitively 
appealing. It is the stretching which generates the relative vorticity. 
The horizontal Laplacian of § transforms directly from cartesian 
to cylindrical coordinates. Equation II-l is therefore the governing 
equation in the latter system also. Let 
E(Zgsr) = R(r)Z (Zo) 


in order to separate variables, and II-l becomes 
f2 


NZ(Zo) 


ZR +ZV2R = 0. TL? 


2. 


For notational convenience, n?(Z.) will be written as No 


Rearranging II-2 yields 
Z V2 
ie. 2 ese ee, II-3 


NG Z, R 


Since a function of Za can be equal to a function of r only if they are 


each constants, then 


_ a) a. 7 
nz : PeBigemsk aa 


Equation II-4 implies 


2 
nz = 1-2 = 0. TiS 
Np 


Equation II-5 must be solved numerically because of the inclusion of 
the Brunt-Vaisala frequency as a function of Zo- Note that the 
selection of sign for the separation constant implies that solutions 
are sought which are exponential in form in the vertical and harmonic 
in the horizontal. Equation II-4 implies 

v2R + k2R = o. IL-6 


Writing the Laplacian in full, one obtains 





R" + R't k*R = 0, [ney 
assuming axisymmetric motion. Equation II-7 is a form of Bessel's 
equation and has the solution 

R =B Jo(kr) + C No (k’). II-8 
No is a Neumann function or Bessel function of the second kind. Since 
the solution must be finite at r = 0, the constant C equals zero. 
With no loss of generality, the radial solution is normalized to unit 
amplitude at r = 0. The separation, then, is the following: 

00 

§ = 2, § (00) Jo (kpr). II-9 
There is an infinity of solutions dependent on the value of k,, the 
arbitrary constant of separation called the wavenumber. In the appli- 
cation of equation II-9 to this study, the solution can be reduced to 
the response to a single, dominant wavenumber. The technique will be 
discussed following a summary of the numerical method used to solve 
equation II-l. 

Equation II-l was solved numerically on the IBM 370/155 computer 
at the University of Miami. A pre-programmed scientific subroutine 
package was utilized which is a modification by Hamming of Milne's 
predictor-corrector method. Since the predictor-corrector method is 
not self-starting, a fourth order Runge-Kutta starting procedure is 
used to generate sufficient points (Ralston, 1965). The numerical 
method requires specification of § and 0€/0Z 6 as initial conditions. 
At the bottom the stretching of an isopycnal is zero. It is assumed 
that a dominant wavenumber for the disturbance exists in the interior. 
It will be shown later that this assumption is realistic. There are, 
then, two arbitrary parameters, wavenumber and slope at the bottom, 


which have to be determined. A solution for specific values of slope 
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and wavenumber can be obtained by devising a scheme which matches the 
numerical solution to the density structure in two dimensions. The 
matching uses up the two degrees of freedom associated with the 
arbitrary parameters. A numerical scheme could be devised which would 
optimize the matching process by some best fit method. It was felt, 
however, that a subjective graphical method was justified considering 
the quality of the data. The approach chosen was to select an inter- 
mediate density surface and vary the input to the computer program 
until the amplitude of the response and the decay of the Jo(kr) Bessel 
function matched the observed perturbation. Essentially, the numerical 
solution, generated as an initial value problem, was forced to meet 
conditions on a surface on which density is conserved. 

The governing equation (II-1) requires the Brunt-Vaisala frequency 
of an undisturbed ocean as input. The input was taken from the 
temperature and salinity tables from R/V ATLANTIS stations 5446 and 
5447 (Fuglister, 1960). The stations were selected because they display 
temperature-salinity characteristics fairly typical of the Sargasso 
Sea and are relatively unperturbed by local variations. Figure 3 shows 


the observed response of the ocean to the storm. 
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Figure 3. Observed mass structure of ocean after passage of Hurricane 
Kara. Isopleths are values of sigma-t (,)- 





IIL. DISCUSSION OF RESULTS 


A. DISCUSSION OF SOLUTION AND 
CONSISTENCY OF RESULTS 

Figure 4 shows a comparison between the computed and observed 
curves. The observed mass distribution has been shifted to a common 
vertical axis. This step is consistent with the assumption of 
axisymmetric motion implied in the transformation of the Laplacian 
in equation II-6 to cylindrical coordinates. As seen in Figure 3, the 
axis of the observed curves slants slightly away from the vertical. 
With the exception of the near surface layer and a few isolated spots 
at the sides, where the assumption of axisymmetric motion is in doubt, 
the difference between the curves is no more than 10% of the total 
perturbation amplitude (trough to crest) of any isopycnal. In most 
cases the percentage is considerably lower. The cross-hatched area 
in Figure 4 indicates the areas where the deviation of the two curves 
is significant. Considerable flattening of the density surfaces occurs 
at least as deep as 230 meters, and there is no evidence of upwelling 
above 100 meters. 

The lack of information about the state of the near surface layer 
prior to and during the storm precludes a complete analysis of the 
causes of the deviation between the two curves. According to Fuglister 
(1947) , average surface temperatures for October in the area of the 
storm range from 25.8°C to 26.1°C. These temperatures correspond very 


well to the surface temperatures found after passage of Kara. If the 
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Figure 4. Comparison between the observed mass structure and the computed 
response after Hurricane Kara. Jo(kr) = 0 at 80 km. 
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near surface layer did not cool significantly, then one must ask if 
upwelled water penetrated the thermocline and was mixed in the upper 
layer. It is possible that the apparent lack of surface cooling in 

the presence of obvious upwelling is the result of either horizontal 
advection in the surface layer moving the cold water away from the place 
of formation,or advection in the interior. The former is more likely, 

but the latter can partially explain the slanting of the axis of 

rotation of the vortex and the fact that the dome is located significantly 
to the south of the published hurricane track. 

Upwelling clearly took place as a response to Kara. One would 
expect to find evidence of downwelling at the sides to satisfy continuity 
requirements. O'Brien (1967) finds downwelling in a numerical model 
and Leipper (1967) finds evidence of it from the changes in shape and 
depth of the thermocline. Surface data taken after Kara do not indicate 
downwelling. The mixed layer depth is remarkably constant across the 
vortex and is typical of the ocean in this region. At the sides of the 
vortex, in the interior, downwelling definitely occurred in balance of 
the upwelling in the center region. This can be seen by the dipping 
of the isopycnals below the undisturbed level. Ths disparity between 
the curves at the top of the pycnocline suggests significant downward 
penetration of momentum below the depth of the mixed layer, leading to 
an effective Ekman depth of about 230 meters. That is, if the Ekman 
layer is considered that region where friction effects dominate, and 
where the conservation equations are not applicable, then the depth 
of significant differences between the observed and computed curves 


may be identified as the Ekman depth. 





It is possible to check the consistency of the results with 
existing ideas of surface wind stress by working backwards from the 
theoretical displacement of the isopycnals. The technique begins with 
the Ekman pumping relation 

Ofw = curl T, III-1l 
where p is density, w is vertical velocity, and T is wind stress. 
Stokes' theorem may be written as 


fr. as = [fcurl T-da, T1i-2 
S A 


or, by substituting equation III-1l for curl T, as 


T-dS = p(lwf-dA. III-3 
pros af 


Assuming an axisymmetric wind stress distribution, equation III-3 
becomes 

owftr2 = 2TrT , III -4 
or 

T = wfr/2. III-5 
If w = 6§/6t where 6§ is the displacement of the most perturbed 
density surface and 5t is the time the storm was stationary, then the 
average stress distribution is 

r= Set III-6 

270t 

The average stress may be found by taking an average over the curve 


g = 50 Jo(kr) . Then 


R 27 
meee 2 | | feed III-7 
2 ee 0 0 2°&t 
es!) i i‘ Jo(kr) -r-do-dr mes 
27 R 00 2° &t 


To simplify the integration, it was decided to utilize the series form 
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of the Bessel function: 


Z 4 6 
a Z + iz Suga, ene 7 
of ) 22(11)2 24(2')2 29 (3!) salad 


Noting that 


and using 
X = kr, dr = dX/k 


to change variables and integration limits, then III-8 becomes 


kr 
_ Eg:f ¢ X Jo(x)dXx 
a an ee III-10 
ae ier 


The dominant wavenumber found by matching the theoretical and observed 
curves is 3.0 x 107°. Therefore, the zero crossing occurs at r = 80 kn. 
Using the well known relation 

T= Pau, = Puvy = CyUTg Pa — 
enables one to use the average stress calculated in equation III-10 to 
find wind velocity, or 


Us = _7T Illicit 
10 Cy:Pa 


where u, is the friction velocity in air, v, is the friction velocity 
in water,P, is the density of air (1.3 x 1073 gm em73), Pp, is the 
density of water (1 em em-3) , Cp is the drag coefficient, and Ujg is 
the wind velocity at a specified height. The wind velocity is highly 
dependent on the choice of drag coefficient. Unfortunately, there are 
numerous estimates of the dependence of drag coefficients on wind speed 
and few agree. From Roll (1965) there seems to be a bunching of 
estimates for near hurricane force winds at Cp = 2.5 x 1073. This is a 


figure O'Brien (1967) also uses and terms "classical."' The results of 


this calculation for various values of the drag coefficient are 
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Eapulated below. Cp = 1.5 x 1073 is a widely quoted value which Malkus 
and Riehl (1960) found worked well in their model. 2.0x 1073 is the 


mean of the two values. 





Cp U10(KTS) 
1.5 x 1073 5) (6 
2.0 x 1073 50.2 
2.5 x 1073 45.2 





Air Force penetrations of the storm during the time it was quasi- 
stationary found wind speeds between 45 and 55 knots. According to 
Colon (1966), the wind distribution in the lower few thousand feet is 
uniform, or nearly so. Therefore the wind speeds should be repre- 
sentative of those near the surface. The tabulated results of the 
calculations are inconclusive. One may say there is no inconsistency 
between the theoretical result and gathered data. Occasionally the 
results of theoretical models are used to fix parameters upon which 
the models critically depend. A definite statement cannot be made in 
this case. However, the results indicate that the drag coefficient 
probably falls within the above range for winds of near hurricane force. 
As will be shown later in the paper, the penetration of the 
response of high wavenumber motion is less than that of lower wave- 
numbers due to the effect of the wavenumber on the exponential-like 
solution in the vertical (equation II-5). The ocean acts as a low- 
pass filter to the response. This means that the difference between 
the two curves at the level where they diverge significantly may be 
due partially to the effect of the response at higher wavenumbers, 


which will change the height of a density surface in some unknown fashion. 
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It is also possible that the presumed center of the vortex was 
not the exact center. If not, then missing the center can account for 
some of the difference between the curves. 
B. THE MIXING PROCESS 
AND MOMENTUM TRANSFER 
Several calculations can be made which place an upper limit on 
the effects of the mixing process. The calculations assume a balance 
in the vertical between the downward diffusion of heat due to mechanical 
mixing and the upwelling of cold water from below the local thermo- 
cline. The governing equation for this process is 
wOT/8Z =k -d27/az2, Teele 
where w is vertical velocity, T is temperature and K« is the eddy 
diffusivity of heat. The solution of III-13 is 
T = Tg exp(w/K)Z. III-14 
K/w is the e-folding distance for the temperature decay and is therefore 
the scale height of the local thermocline in this problem. Hy, the 
thermocline scale height, may be found by taking the logarithm of III-14, 
obtaining the following: 
Of = UnIp + Z/Hr. III-15 
Plotting #f€ vs Z and finding the slope determines the scale height. 
Equation III-15 is the equation of a straight line on semi-log paper, 
and 1/Hp is the slope of the line. Figure 5 is such a plot. Since 
pure mixing is assumed, the vertical velocity may be approximated by 
5h/&t, the difference in lifting between the observed and calculated 
curves divided by the relevant time scale of one day. k« determined in 


this manner is 1400 em /sec. This value of K is several orders of 
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magnitude too large for abyssal processes but is probably the correct 
order for surface processes where eddy viscosity and eddy diffusivity 
may be nearly the same. If the point of divergence of the curves is 
at 230 meters, the eddy viscosity found from the Ekman depth formula, 
dg = mV 2V/£ » is 2000 cm/sec. 1400 cm2/sec is almost certainly an 
overestimate of the eddy diffusivity since Ekman divergence effects 
are neglected. If the depth of the mixed layer is neglected in the 
calculation of the eddy viscosity, dg = 180 meters and v = 1200 cm/sec. 
This calculation assumes the existence of two separate regimes: the 
mixed layer where eddy viscosity coefficients are very high; and the 
layer of Ekman divergence below the frictional layer, indicated by the 
disparity between observed and computed results. The similarity between 
the values of eddy viscosity and eddy diffusivity is presented for 
completeness and to suggest topics for further investigation. 
Under conditions of negligible advection, it is possible to 

extend an energy calculation by Turner (1969). Assuming pure mixing, 
the change in potential energy between the observed and theoretical 
curves can be used to work backwards to find wind velocity. Selecting 
an isopycnal where the two curves do not diverge and calling the density 
Po a reference density, one can plot Pz7-fPqg vs z2 (depth squared). 
Taking one half the area under the resulting curve solves the energy 
equation, 

Epot = 8) (Pz-Po)Z ad Z, TII-16 
graphically. Figure 6 is a plot of 4p vs z* for the center station, 
where maximum upwelling occurred. The change in potential energy found 


planimetrically is 10.14 x 10/ ergs cm72. By dimensional arguments, 
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Figure 6. Density difference (zm/cm?) plotted vs square of depth (cm?) 
for potential energy calculation. Station #24 is compared to the center 
of the computed curves. 
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Turner asserts that it is reasonable to write the time rate of change 


of energy as 


—For meee utD = C>pU19°. TII-17 
Ct 


The tabulation below gives values of wind speed for two values of the 


constant Co found by Turner. 


Cy U1 (KTS) 
2.6 x 1078 7b 
1 ator” 90 


Both values of wind speed are too high if pure mixing is responsible 

for the potential energy change. If one accepts 50 knots as a 
reasonable value of wind speed during the storm, then surface mixing 
accounts for about half the energy change. That is, at least half of 
the density field discrepancy from the theoretical curves must derive 
from the deep Ekman divergence. As the drag coefficient, Turner's 
constants are not well defined. The dependence of Turner's constants 
and the drag coefficient on wind speed is unknown. Turner's method 
requires a one-dimensional process and a regime in which surface cooling 
and resultant convective overturning are unimportant. Both processes 
are usually significant under hurricane conditions, although their 
effect may be minimal in this case since surface cooling apparently did 
not occur. Selecting the center station minimizes the advective effects, 
since § is a maximum there and 90§/6ér approaches zero. As Turner 


correctly points out, his method is illustrative but hardly conclusive. 


C. THE EFFECT OF LENGTH SCALE ON RESPONSE 


Some interesting properties of the solution as a function of scale 
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of forcing and depth are illustrated in Figure 7. The governing equation 
(II-1) was solved using the same initial conditions and density profile 
but different values of wavenumber. The length scale was determined 
by selecting the value of radius which corresponded to the first zero 
of the Jo Bessel function, found from the formula (1/wavenumber) x 2.4. 
The values of § were normalized to 100 meters of Ekman pumping at 
Zo = 25 meters and isopleths of constant E were plotted. The response 
may be split into a barotropic and baroclinic regime. In the barotropic 
regime, response is characterized by the stretching being a linear 
function of depth, since stretching is proportional to vertical velocity. 
The transition from a baroclinic to a barotropic response is rather 
broad in scale, occurring between 300 to 500 km. For scales greater 
than 500 km, the response may be continued to be barotropic. The 
response for high wave numbers is confined to a shallow surface layer. 
The latter provides the rationale for basing the solution of the 
vorticity equation on the existence of a dominant wavenumber in the 
interior. The dependence of the vertical structure of the response on 
the scale of the motion suggests a definite problem in the representation 
of a continuously stratified ocean by a two-layer model. For the higher 
wavenumbers, the two layer model lacks sufficient resolution to describe 
the response adequately. 
The ordinate in Figure 7 is undisturbed depth (Zo). Since depth 
is equal to the sum of the undisturbed depth and the stretching, zero 
depth in the figure should not be confused with the free surface. The 
free surface lies considerably below Zo = 0. 


It is anticipated that similar diagrams can be generated using 
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Figure 7. Isopleths of stretching, E(m), as a function of depth and 
length scale. Values of § are normalized to 100 meters pumping at 
25 meters matching depth. 
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density distributions typical of other areas of the ocean. The diagrams 
should prove useful in predicting the response in these areas. A 
density profile with higher resolution than that obtained from 


Fuglister's Atlantic Ocean Atlas is preferable. 
D. GEOSTROPHIC VELOCITIES 


It is possible to calculate geostrophic velocities using the 
solution of equation II-l. Noting 
Z(p or) = ZQ(p) oP §(Z 94) 5) ITI-20 


and from Appendix A, 


Qa Bz) Bp 21 
8x) Z0 | Zo ” nee 


one may use the thermal wind equation, 


OU = -g op , igigeryy 
0Zo of drj/Zo 


to find velocities throughout the column and for any radius. 


Substituting ILII-21 into III-22, one obtains 


wu _g Z\ op . LIL -23 
Zo $f z), Zo 
Equation III-23 may be integrated in the vertical: 
Ur} Zc 
0 
in du =i & oe . (OZ dZ -). III-24 
Uz 5 f 7" — OZ9 Or /P 
Noting Uz (bottom) = 0 and No = =e Op , equation III-24 becomes 
! 
Z0 
1 2 [eZ 
Oe a N dZp . TIT-25 
“0 ae s (2), £ 


Introduce III-20 for Z in III-25, obtaining 
20 
» 1 2 ,(wOk 
Uzi = -= Na *\= aan. ITI-26 
Zo f !. 0 (28) 0 


Geostrophic velocities were calculated numerically using equation III-26 
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and the result established previously, § =  Jo(kr). The derivative 
Of = “is 

Seer = lear auia(kx) TII-27 
Equation III-26 is a valid expression for the velocities induced by 
Kara, except in the Ekman layer, where the predicted density field 
departs from the observed. Veronis (1956) showed that motion which 
has forcing on scales larger than the baroclinic radius of deformation 
is essentially geostrophic. In this case, the scale of forcing, 80 km, 
is more than twice the baroclinic radius of deformation, 36 km. 
Figure 8 is a comparison of the geostrophic velocities based on 
equation III-26 and velocities calculated using the dynamic method 
(called the observed velocities). The reference level for calculating 
the observed velocities was found using Defant's method. The level 
was conspicuously present in all the profiles. Although the percentage 
of error is high at the lower depths, the absolute magnitude difference 
is low. In general, the magnitude of error remains the same as depth 
decreases. The two velocity profiles agree rather well in magnitude 
and general shape. Aside from the reasons for the difference between 
the curves in the Ekman layer, which have been discussed, the calculation 
of geostrophic velocities poses special problems. In this case, the 
distances between stations used to calculate the velocities by the 
dynamic method were extreme. It was found that velocities at a given 
radius and depth varied as much as 5 cm/sec depending on the particular 
combination of stations used. The curves represent the smoothest 
combination. The observed velocities were computed from the non-centered 
observed data, which certainly accounts for some of the difference 


between the curves. The effect of centripetal acceleration has been 
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Figure 8. Comparison of numerically calculated geostrophic velocities 
(cm/sec) with velocities found using the dynamic method (observed). 
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neglected in the calculation of velocities, and in the formulation of 
the vorticity conservation model. The inclusion of centripetal 
accelerations makes the model non-linear, introducing an undesired 
degree of difficulty in solving the vorticity equation. A measure of 


the importance of the centrifugal force term is the Rossby number, _Y 


(Neumann and Pierson, 1966). Choosing r = 40 km, v = 40 cm/sec, and 
f =7.5 x 107 sec7l, the Rossby number is 0.13. It can be shown that 
geostrophic velocity is an overestimate of the gradient velocity if 
curvature is important (Haltiner and Martin, 1957). 

The success in obtaining substantial agreement between the 
theoretical mass field and velocities and the observations supports the 


validity of the initial model formulation. 





. IV. SUMMARY AND CONCLUSIONS 


A model based on the conservation of potential vorticity was 
solved numerically. The solution was utilized to investigate some 
aspects of the oceanic response to a tropical storm as a function of 
horizontal scale. In particular, the interior response of a moderately 
intense, quasi-stationary cyclone was duplicated. The matching was 
successful below a depth of 230 meters which was associated with the 
Ekman depth. The lack of information concerning the state of the 
ocean in the near-surface layer prior to the passage of Hurricane Kara 
compounded the problem of discovering what caused the divergence of 
the observed and computed curves. The divergence was attributed to the 
neglect of high wavenumber modes when solving for the interior response , 
and to some ill-defined combination of vertical mixing and momentum 
transfer by friction. Gross estimates indicate that the two processes 
are of equal importance in this particular case. An upper limit for 
the eddy diffusivity is 1400 cm*/sec. Vigorous mixing, if it occurred, 
may not have resulted in significant cooling of the mixed layer. 

Geostrophic velocities were calculated from the raw data using 
Defant's method to find a reference level of no motion. These 
velocities were compared to those calculated from the numerical solution. 
The velocities agree rather well below the near-surface layer. Since 
the velocity computation integrates the stretching, the comparison of 
velocities should not be expected to be as successful as the comparison 
of € itself. The integration of § means that errors are accumulated 
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over the entire column. Geostrophic vertical shear matches better 
than geostrophic velocities. It is plausible that the nearly constant 
error in magnitude between the velocity fields is a result of 
inaccuracies in the dynamic computations; specifically, limitations in 
the Defant method of obtaining the level of no motion may be the cause. 

It was found that the scale of forcing can have dramatic effects 
on the penetration of the response. For equal Ekman pumping, higher 
wavenumber (short scale) response does not penetrate very deep into 
the ocean compared to low wavenumber response. Thus the scale of the 
disturbance must be considered when choosing the depth and number of 
layers for modelling a continuously stratified fluid. A two layer 
model may lack the resolution necessary to accurately describe the 
response in the upper layer. Adequate modelling of the near-surface 
layer by a model such as the one uSed in this study requires the 


inclusion of high wavenumber motion. 


APPENDIX A 


THE NATURAL COORDINATE SYSTEM 


si 
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Although the standard cartesian coordinate system is the one 
normally used in hydrodynamics, there is no reason why the vertical 
coordinate must be geometric height, a fact recognized long ago by 
meteorologists and some oceanographers. The use of a "natural" vertical 
coordinate, such as pressure, density or potential density, can simplify 
an analysis but render boundary conditions formidable. 

In the potential vorticity conservation model, density will be 
utilized as the "natural" vertical coordinate. This is advantageous 
because it enables one to solve directly for the displacement of an 
isopycnal from equilibrium. The method of coordinate transformation is 
discussed in Haltiner and Martin (1957). It is summarized below since 
the method is implied in the model and is used directly in deriving the 
geostrophic velocities (equation III-26). 

Consider $(x,y,z,t) = constant. The change of 1\(x,y,z,t) along 
® = constant is equal to the change along the horizontal plus the 


change along the vertical due to the slope of ¢, i.e. From Figure 9, 


om) - iM) id 1) ee aa 
Now, 
= ON\ 4 + SN) 5 A-2 
2) Ox ¥ 52st oz neo 


by the chain rule at a constant time, tg. A-2 may be written as 


oT = 31) + 2) Oz! ; A-3 
Ox y 9,tg Ox y»Z2,to Oz X»Y,to Ox y9,to 


where Oz is the slope of the ¢ surface in the x direction. Let Tl] =p 
x 


and A-3 becomes 


OQ 
N 


Op\ . 9p) _ op . oz 
ox /z Ox/G oz Ox 


: 
/ 
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Figure 9. % = constant surface separated into horizontal and vertical 
components to illustrate coordinate transformation technique. 





APPENDIX B 


THE POTENTIAL VORTICITY CONSERVATION MODEL 
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The governing equation for the interior motion of a hurricane, it 
is hypothesized, expresses the conservation of potential vorticity. 
The equations of motion used are non-linear, time dependent, hydrostatic 
equations in a rotating system characterized by the Coriolis parameter, 
£, which is a function of the horizontal coordinates alone: 


0 B-1 


Up a tleme Vu - Exe + We 


caterer vy + fu + we 0; B-2 


e 


where T = . + eZ. The vorticity equation may be formed by cross- 


differentiating equations B-l and B-2, resulting in B-3: 


at Ox 


B-3 may be written as 


+ u Btvd (vy-uy) tft, -uy) (uxt) + uf, + vfy = 0 B-3 


eC) 3 3 = 
ce + ty -u_)+ : +v.) =0, B-4 
E ad | (itv uy) (iv,, uy) Cu, w 
since f£f is not a function of time. That is, 
d = 
aq = an _ = 0, B-5 
ae Civ uy) Cine. uy) (uty) 0 
where 
d 3 re) 3 
—— aa + a . 
dé wer “ox | By P 
The divergence equation may be formed by taking 5, (equation B-4) 


and adding to 8. (equation B-2). 


2 yey u tuv.. = 0 B-6 


d yen a wn 
a = -— - 
dt (uy y) f(v., uy) fyu-fy, Sy Ks oy x 


y 
Assuming geostrophic motion reduces the balance equation to 
f(v,-u_) = V2n. Ba 
(We ” 
The continuity equation is 
= -8 
Sez 0. B 


In the vertical, w depends only on z(p), therefore 
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Cu, tv a + Ww. B=, 
dp az 


Noting w = dz/dt, equation B-8 becomes 


Se uytyy) +¢ 0z\ = 9, B-9 
t\ dp, 
substitute equation B-9 into the vorticity equation (B-4) to obtain 
d dz /a 
dé iy mam) =e(f£ Fv, - u.) - deer B-10 
dt x y x y dz /dp 
Rearranging terms, equation B-10 becomes 
SG + vx 7 Uy) S<(8z /8p) 7 3-11 
See Ne = a = 3 - 
P+ vy - uy Oz / Ap 
or, by definition of the logarithm, 
Vem U. tt 
g. (eet) 0. B-12 
6z/ap 


Equation B-12 implies that for each fluid element 


Vy Tuy tt = constant. B-13 
Oz /dp 


The existence of a uniform reference state implies that the constant is 
a function of density alone. If the reference state is rest, relative 
vorticity is zero, and the constant, $(p), equals £/(8z/8p) 9. The 


vorticity equation, then, becomes 


[ctl . -£ . B-14 


Substituting equation B-7 into B-l2 and rearranging, one obtains 
(HE 2)(22) = oz Bos 

£2 /\ap/0 ap 
Next, apply the hydrostatic equation, 

=PB4 2Z. B-16 

p 

Rearranging and taking the partial derivative with respect to z results 
in 


3 
5, (T-gz)p = PL. B17 





bh) 


Assuming hydrostatic equilibrium, B-17 becomes 


2 (m-gz)p = -gp. B-18 
oz 
Integrating B-18 results in the following: 
Z 
pTI-Ppgz = Pg To - i ogdz. B-19 
0 
The last term in equation B-19 can be integrated by parts: 
Z p 
-{ogdz = -pgz + |gzdp - B-20 
0 Po 


Substituting B-20 into B-19, the desired forms of the hydrostatic 


equation is achieved: 


p 
T= To += Jezdp B-21 
PO 
Take the Laplacian of B-21: 
ee 2h fae? 
6 + er eee, B22 
0 


Substitute equation B-22 into equation B-15: 


p 
ez Oz = l (v2n + l v2 d (#) 
dp (=), £2 ‘ Js 8p/o 


Therefore, if Pop is chosen as the surface density, where h is the free 


surface height, 


p 
oz Bs oz = E—(92h + L reap) (32) B-23 
dp = 0 f2 9 J 3p 0 


The depth, z, may be split into two parts, 


Z = ifp) + E(x sy »Z 0) ° 


Equation B-23 becomes, then, after rearranging, 
p 
2 ee. = £2 ae 4 =9h 
fe) se g op (3 0 


Next, take the derivative of B-24 with respect to p: 


2 
lo2e = £°po /9§ 3°] - f2p99 0 /o9§ B-25 
=—V/ ee ° 
p : 890 (38) ($8 0 gp0 oP Oz 
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Recall 
a _ (a2) a 


so equation B-25 becomes 


1 ye = -f4 a¢E 
x ¥ ; 
PO : pONZ(Zo) aZo2 


or 


ve +2 2B 


ay ° B-26 


Equation B-26 is the governing equation used in the analysis described 
in the text. The basic derivation is attributed to Rooth (1970). 


Several steps were expanded by this author for clarity and completeness. 
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